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ABSTRACT 

We  report  on  the  analysis  of  the  first  complete  far-infrared  spectrum  (43-197  pm)  of  the  Seyfert  2  galaxy 
NGC  1068  as  observed  with  the  Long  Wavelength  Spectrometer  (LWS)  on  board  the  Infrared  Space  Observatory 
(ISO).  In  addition  to  the  seven  expected  ionic  fine-structure  emission  lines,  the  OH  rotational  lines  at  79,  119,  and 
163  pm  were  all  detected  in  emission,  which  is  unique  among  galaxies  with  full  LWS  spectra,  where  the  119  pm 
line,  when  detected,  is  always  in  absorption.  The  observed  line  intensities  were  modeled  together  with  ISO  Short 
Wavelength  Spectrometer  (SWS)  and  optical  and  ultraviolet  line  intensities  from  the  literature,  considering  two 
independent  emission  components:  the  active  galactic  nucleus  (AGN)  component  and  the  starburst  component 
in  the  circumnuclear  ring  of  ~3  kpc  in  size.  Using  the  UV  to  mid-lR  emission  line  spectrum  to  constrain  the  nu¬ 
clear  ionizing  continuum,  we  have  confirmed  previous  results:  a  canonical  power-law  ionizing  spectrum  is  a 
poorer  fit  than  one  with  a  deep  absorption  trough,  while  the  presence  of  a  “  big  blue  bump”  is  ruled  out.  Based  on 
the  instantaneous  starburst  age  of  5  Myr  constrained  by  the  Br7  equivalent  width  in  the  starburst  ring,  and  star- 
burst  synthesis  models  of  the  mid-  and  far-infrared  fine-structure  line  emission,  a  low-ionization  parameter  (U  = 
1CL3-5)  and  low  densities  ( n  =  100  cnr3)  are  derived.  Combining  the  AGN  and  starburst  components,  we  suc¬ 
ceeded  in  modeling  the  overall  UV  to  far-lR  atomic  spectrum  of  NGC  1068,  reproducing  the  line  fluxes  to  within 
a  factor  of  2.0  on  average  with  a  standard  deviation  of  1 .3,  and  the  overall  continuum  as  the  sum  of  the  contribu¬ 
tion  of  the  thermal  dust  emission  in  the  ionized  and  neutral  components.  The  OH  119  pm  emission  indicates 
that  the  line  is  collisionally  excited  and  arises  in  a  warm  and  dense  region.  The  OH  emission  has  been  modeled 
using  spherically  symmetric,  nonlocal,  non-LTE  radiative  transfer  models.  The  models  indicate  that  the  bulk 
of  the  emission  arises  from  the  nuclear  region,  although  some  extended  contribution  from  the  starburst  is  not  ruled 
out.  The  OH  abundance  in  the  nuclear  region  is  expected  to  be  ~  1 0  5 ,  characteristic  of  X-ray-dominated  regions. 
Subject  headings:  galaxies:  active  —  galaxies:  individual  (NGC  1068)  —  galaxies:  nuclei  — 
galaxies:  Seyfert  —  galaxies:  starburst  —  infrared:  galaxies 


1.  INTRODUCTION 

NGC  1068  is  known  as  the  archetypical  Seyfert  type  2  gal¬ 
axy.  It  is  nearby  and  luminous  (Gr  =  2  x  1 0 1 1  Le ;  Bland-Hawthom 
et  al.  1997),  and  it  has  been  extensively  observed  and  studied  in 
detail  from  X-rays  to  radio  wavelengths.  With  a  measured  red- 
shift  ofz  =  0.0038  (Huchra  et  al.  1999;  corresponding  to  a  dis¬ 
tance  ofZ>=  15.2  Mpc  for  Ho  =  75kms~1  Mpc~*),  itprovides 
a  scale  of  only  ~74  pc  arcsec-1.  A  central  nuclear  star  cluster 
has  an  extent  of  (Thatte  et  al.  1997),  and  a  2.3  kpc  stellar 
bar  observed  in  the  near-IR  (Scoville  et  al.  1988;  Thronson  et  al. 
1989)  is  surrounded  by  a  circumnuclear  starburst  ring.  Telesco 

1  ISO  is  an  ESA  project  with  instruments  funded  by  ESA  Member  States 
(especially  the  PI  countries:  France,  Germany,  the  Netherlands,  and  the  United 
Kingdom)  with  the  participation  of  ISAS  and  NASA. 


et  al.  (1984)  found  that  the  infrared  emission  in  NGC  1068  was 
due  to  both  the  Seyfert  nucleus  (which  dominates  the  10  pm 
emission)  and  to  the  star-forming  regions  in  the  bright  ^3kpc 
circumnuclear  ring  (which  emits  most  of  the  luminosity  at  X  > 
30  pm).  A  Bry  imaging  study  (Davies  et  al.  1998)  showed  a 
similar  morphology  and  indicated  that  a  short  burst  of  star  for¬ 
mation  occurred  throughout  the  circumnuclear  ring  of  1 5"- 1 6" 
in  radius  within  the  last  4-40  Myr.  CO  interferometer  observa¬ 
tions  revealed  molecular  gas  very  close  to  the  nucleus  (~0"2), 
suggesting  the  presence  of  ~  10s  M0  within  the  central  25  pc 
(Schinnerer  et  al.  2000).  Recent  high-resolution  H2  line  emis¬ 
sion  mapping  indicates  the  presence  of  two  main  nuclear  emis¬ 
sion  knots  with  a  velocity  difference  of  140  km  s_1,  which,  if 
interpreted  as  quasi-Keplerian,  would  imply  a  central  enclosed 
mass  of  10s  M0  (Alloin  et  al.  2001). 
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In  this  paper  we  present  the  first  complete  far-infrared  spec¬ 
trum  from  43  to  197  pm  showing  both  atomic  and  molecular 
emission  lines  (§  2).  We  model  the  composite  UV  to  far-lR 
atomic  emission-line  and  continuum  spectrum  from  our  data  and 
the  literature,  using  photoionization  models  of  both  the  active 
nucleus  and  the  starburst  component  (§  3).  We  also  model  the 
mid-  to  far-lR  continuum  emission  using  a  radiative  transfer  code 
and  graybody  functions  for  the  neutral  molecular  components 
(§  4).  Moreover,  two  different  nonlocal,  non-LTE  radiative  trans¬ 
fer  codes  have  been  used  to  model  the  OH  lines  (§  5).  Our  con¬ 
clusions  are  given  in  §  6. 

2.  OBSERVATIONS 

NGC  1068  was  observed  with  the  Long  Wavelength  Spec¬ 
trometer  (LWS;  Clegg  et  al.  1996)  on  board  the  Infrared  Space 
Observatory  ( ISO ;  Kessler  et  al.  1996),  as  part  of  the  Guaran¬ 
teed  Time  Programme  of  the  LWS  instrument  team.  The  full  low- 
resolution  spectrum  (43-197  pm)  of  NGC  1068  was  collected 
during  orbit  605  (1997  July  13).  Two  on-source  full  scans  (15,730  s 
of  total  integration  time)  and  two  off-source  (6'  north)  scans 
of  the  [C  n]  158  pm  line  (3390  s  of  total  integration  time)  were 
obtained.  On-  and  off-source  scans  had  the  same  integration 
time  per  spectral  step.  Because  of  the  design  of  the  LWS,  a  short 
spectral  scan  of  equal  sensitivity  to  the  on-source  spectrum  was 
obtained  simultaneously  with  the  158  pm  data  at  sparsely 
spaced  wavelengths  across  the  LWS  range. 

The  LWS  beam  is  roughly  independent  of  wavelength  and 
equal  to  about  80".  The  spectra  were  calibrated  using  Uranus, 
resulting  in  an  absolute  accuracy  better  than  30%  (Swinyard 
et  al.  1996).  The  data  analysis  has  been  done  with  the  ISO  Spec¬ 
tral  Analysis  Package  (ISAP),2  starting  from  the  autoanaly¬ 
sis  results  processed  through  the  LWS  version  7-8  pipeline 
(1998  July).  To  be  confident  that  newer  versions  of  the  pipeline 
and  calibration  files  did  not  yield  different  results,  we  have 
compared  our  data  with  the  results  obtained  using  pipeline  10.1 
(2001  November),  and  we  did  not  find  significant  differences  in 
the  line  fluxes  or  the  continuum. 

All  the  full  grating  scans  taken  on  the  on-source  position  and 
the  two  sets  of  data  on  the  off-source  position  were  separately 
co-added.  No  signal  was  detected  in  the  off-source  co-adds.  The 
emission-line  fluxes  were  measured  with  ISAP,  which  fits  poly¬ 
nomials  to  the  local  continuum  and  Gaussian  profiles  to  the 
lines.  In  all  cases  the  observed  line  widths  were  consistent  with 
the  instrumental  resolution  of  the  grating,  which  was  typically 
1500  km  s_1.  The  integrated  line  fluxes  measured  independently 
from  data  taken  in  the  two  scan  directions  agreed  very  well,  to 
within  10%.  The  on-source  LWS  spectrum  that  resulted  from 
stitching  the  10  LWS  channels  together  using  small  multiplica¬ 
tive  corrections  in  order  to  match  the  overlapping  regions  of  each 
channel  with  its  neighbors  is  shown  in  Figure  4.  LWS  spectra 
of  sources  that  are  very  extended  within  the  instrument  beam 
or  that  peak  off-center  are  typically  affected  by  channel  fringing 
in  the  continuum  baseline  (Swinyard  et  al.  1998).  Fortunately, 
these  spurious  ripples  are  hardly  noticeable  in  our  LWS  spectrum, 
presumably  because  the  far-IR  continuum  is  centrally  concen¬ 
trated  toward  the  center  of  the  LWS  80"  beam. 


2  The  ISO  Spectral  Analysis  Package  is  a  joint  development  by  the  LWS  and 
SWS  Instrument  Teams  and  Data  Centers.  Contributing  institutes  are  Centre 
d’Etude  Spatiale  des  Rayonnements  (France),  Institute  d’Astrophysique  Spatiale 
(France),  Infrared  Processing  and  Analysis  Center  (United  States),  Max-Planck- 
Insitut  fur  extraterrestrische  Physik  (Germany),  Rutherford  Appleton  Laboratories 
(United  Kingdom),  and  the  Space  Research  Organization  (Netherlands). 


Besides  the  LWS  observations,  we  also  use  the  Short  Wave¬ 
length  Spectrometer  (SWS)  observations  presented  by  Lutz  et  al. 
(2000),  to  extend  the  wavelength  and  ionization-level  coverage. 
Table  1  presents  all  the  ISO  line  flux  measurements  including 
those  from  the  SWS  with  their  respective  aperture  sizes. 

3.  THE  FINE  STRUCTURE  LINES 

To  be  able  to  better  constrain  the  modeling  of  the  line  emission 
of  NGC  1068,  we  have  combined  our  far-infrared  fine- structure 
line  measurements  (Table  1)  with  ultraviolet,  optical,  and  infrared 
spectroscopic  data  from  the  literature  (Kriss  et  al.  1992;  Marconi 
et  al.  1996;  Thompson  1996;  Lutz  et  al.  2000).  The  complete 
emission-line  spectrum  of  NGC  1068  from  the  ultraviolet  to  the 
far-IR  includes  several  low-ionization  lines  that  are  primarily 
produced  outside  the  narrow-line  region  (NLR)  of  the  active 
nucleus,  as  well  as  intermediate-ionization  lines  that  originate 
from  both  starburst  and  active  galactic  nucleus  (AGN)  emis¬ 
sion.  For  this  reason,  we  find  that  no  single  model  satisfactorily 
explains  all  the  observed  emission  lines.  We  identify  two  main 
components: 

1.  An  AGN  component  (the  NLR),  exciting  the  high- 
ionization  lines  and  contributing  little  to  the  low-  to  intermediate- 
ionization  lines. 

2.  A  starburst  component  in  the  circumnuclear  ring  of  the 
galaxy  (e.g.,  Davies  et  al.  1998)  that  produces  the  low-ionization 
and  neutral  forbidden  lines  and  some  of  the  emission  in  the 
intermediate-ionization  lines.  This  component  should  also  pro¬ 
duce  emission  associated  with  photodissociation  regions  (PDRs; 
e.g.,  Kaufman  et  al.  1999)  at  the  interface  with  the  interstellar 
medium  of  the  galaxy. 

In  this  section  we  examine  separately  the  two  components 
that  produce  the  total  fine-structure  emission-line  spectrum  of 
NGC  1068,  namely,  the  AGN  and  the  starburst,  for  which  we 
propose  two  different  computations,  and  we  add  together  these 
components  to  reproduce  the  overall  observed  spectrum  from 
the  UV  to  the  far-IR  in  §  3.3. 

3.1.  Modeling  the  AGN 

The  first  photoionzation  model  predictions  of  the  mid-  to  far- 
infrared  emission-line  spectra  of  the  NLRs  of  active  galaxies 
were  presented  by  Spinoglio  &  Malkan  (1992),  well  before  the 
ISO  observations  could  be  collected.  Alexander  et  al.  (2000) 
used  the  observed  high-ionization  emission  lines  to  model  the 
obscured  ionizing  AGN  continuum  of  NGC  1068  and  found  that 
the  best-fit  spectral  energy  distribution  (SED)  has  a  deep  trough 
at  4  ryd,  which  is  consistent  with  an  intrinsic  “big  blue  bump” 
that  is  partially  obscured  by  ~6  x  1019  cm~2  of  neutral  hydrogen 
interior  to  the  NLR.  Following  their  results,  we  have  simulated 
their  models,  although  using  a  different  photoionization  code, 
CLOUDY  (ver.  96;  Ferland  et  al.  1998;  Ferland  2000),  and  we 
have  varied  the  shape  of  the  ionizing  continuum  to  include  the 
ionizing  continuum  derived  in  Pier  et  al.  (1994).  Our  goal  was  to 
test  whether  the  Alexander  et  al.  (2000)  results  were  unique  and 
to  fit  the  remaining  emission  by  a  starburst  component,  thereby 
deriving  a  composite  model  of  the  complete  emission-line  spec¬ 
trum  of  NGC  1068. 

Specifically,  we  explore  three  plausible  AGN  SEDs.  Model 
A  assumes  the  best-fit  ionizing  spectrum  derived  by  Alexander 
et  al.  (2000),  i.e.,  with  a  deep  trough  at  4  ryd  (log/  =  -27.4, 
—29.0,  —27.4,  and  —28.2  at  2,  4,  8,  and  16  ryd,  respectively). 
An  intrinsic  nuclear  spectrum  of  NGC  1068  has  also  been 
inferred  by  Pier  et  al.  (1994).  Model  B  assumes  the  original 
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TABLE  1 


Measured  Line  Fluxes 

FROM  THE  LWS  AND 

SWS  Grating  Spectra,  with  1 

cr  Uncertainties 

X 

Flux 

Aperture 

Line 

Otm) 

(10-13  ergs  s-1  cm-2) 

(arcsec2) 

Reference 

[Si  i\  7’..  •  71.  . 

2.584 

3.0 

14  x  20 

i 

[Mg  vm]  2Py 2  ->  2P1/2 . 

3.028 

11  ±  1.1 

14  x  20 

i 

[Si  IX]  7'  >  7',:  . 

3.936 

5.0  ±  0.6 

14  x  20 

i 

[Mg  iv]  2P1/2  ->  2Py 2 . 

4.487 

7.6  ±  1.5 

14  x  20 

i 

[Ar  vi]  2P]/2  2P\ /2 . 

4.529 

15  ±  3 

14  x  20 

i 

[Fe  n]  a4Fi)/2  -  a6D9/2 . 

5.340 

5.0 

14  x  20 

i 

[Mg  vii]  7’2  •  7’; . 

5.503 

13 

14  x  20 

i 

[Mg  v]  7>.  7>2 . 

5.610 

18  ±  2 

14  x  20 

i 

[Ar  ii]  2Pl/2  ->  2Py 2 . 

6.985 

13 

14  x  20 

i 

[Na  ill]  2Pl/2  2Pi/2 . 

7.318 

5.8 

14  x  20 

i 

[Ne  vi]  2P3/2  ->  2 Pi/ 2 . 

7.652 

110  ±  11 

14  x  20 

i 

[Fe  vii]  3 F4  — >  3F3 . 

7.815 

3.0 

14  x  20 

i 

[Ar  v]  7’..  •  7’, . 

7.902 

<12 

14  x  20 

i 

[Na  vi]  3P '2  ->  3Pi  . 

8.611 

<16 

14  x  20 

i 

[Ar  hi]  3P,  ->  3  Pi . 

8.991 

23.0  ±  3.3 

14  x  20 

i 

[Fe  vii]  3F3  — >  3 Ft . 

9.527 

4.0 

14  x  20 

i 

[S  iv]  2P3/2  -t-  2P i/2 . 

10.510 

58  ±  6 

14  x  20 

i 

[Ne  n]  2P2/2  — >  2P i/2 . 

12.813 

70 

14  x  27 

i 

Ar  v]  7’  -7’,. . 

13.102 

<16 

14  x  27 

i 

[Ne  v]  ->  3 P{ . 

14.322 

97  ±  9.7 

14  x  27 

i 

[Ne  iii]  3P\  -»  3  P-, . 

15.555 

160  ±  32 

14  x  27 

i 

[Fe  n]  “4F7/ 2  -  “4F9/2 . 

17.936 

<10 

14  x  27 

i 

Sh7’2  -7>  . 

18.713 

40 

14  x  27 

i 

[Ne  v]  7'  ->  7’:, . 

24.317 

70  ±  7 

14  x  27 

i 

[O  iv]  2p,/ 2  -»  2P1/2 . 

25.890 

190  ±  20 

14  x  27 

i 

[Fe  n]  “7)  2  -  a6D9  2 . 

25.988 

8 

14  x  27 

i 

;s  m  7’,  -► 3 Po . 

33.481 

55 

20  x  33 

i 

[Si  n]  2Py 2  ->  2P\/2 . 

34.814 

91 

20  x  33 

i 

[Ne  iii]  3 P0  ->  7’i  . 

36.013 

18 

20  x  33 

i 

[O  iii]  7’-  ->•  7'  . 

51.814 

114  ±  3 

80 

2 

[N  III]  2Py 2  ->  2 P\ /2 . 

57.317 

51.4  ±  2.5 

80 

2 

[O  i]  7',  ->  3P2 . 

63.184 

156  ±  1 

80 

2 

[O  III]  7’  >  7’:,  . 

88.356 

111  ±  1 

80 

2 

[N  ii]  3 p2  ->  3  P\ . 

121.897 

30.5  ±1.1 

80 

2 

[O  1]  7c  -►  7’  . 

145.525 

11.9  ±  0.4 

80 

2 

[c  „]  2PV2  -I  2p1/2 . 

157.741 

216  ±  1 

80 

2 

OH  2n]/25/2-2n3/,3/2 . 

34.60/34.63 

<3 

20  x  33 

2 

OH  2n1/2l/2-2II1/,3/2 . 

79.11/79.18 

14.4  ±  1.5 

80 

2 

OH  2n3/25/2-2lI3/23/2 . 

119.23/119.44 

11.9  ±  1.2 

80 

2 

oh  2n1/23/2-2nl/2i/2 . 

163.12/163.40 

7.42  ±  0.65 

80 

2 

References. — (1)  From  Lutz  et  al.  2000  and,  where  errors  are  available,  Alexander  et  al.  2000;  (2)  this  work. 


ionizing  spectrum  derived  from  Pier  et  al.  (1994).  Model  C 
assumes  an  SED  with  a  big  blue  bump  superposed  on  the  Pier 
etal.  (1994)  ionizing  continuum  (log/  =  —25.8,-25.8,-25.8, 
and  —27.4  at  2,  4,  8,  and  16  ryd,  respectively),  as  expected  for 
the  thermal  emission  of  an  accretion  disk  around  a  central  black 
hole.  These  three  AGN  ionizing  continua  are  plotted  in  Fig¬ 
ure  1.  For  each  of  the  models  A,  B,  and  C,  we  have  used  two 
component  models  with  the  same  parameters  as  in  Alexander 
et  al.  (2000):  component  1  has  a  constant  hydrogen  density  of 
104  cnr3,  an  ionization  parameter  U  =  0.1,  a  covering  factor 
c  =  0.45,  a  filling  factor  of  6.5  x  10~3  with  a  radial  dependence 
of  the  form  r~2,  and  extends  from  ~21  to  ~109  pc  from  the 
center;  component  2  has  a  density  of  2xl03  cm'3,  an  ioniza¬ 
tion  parameter  U  =  0.01,  a  covering  factor  of  c  =  0.29,  a  fill¬ 
ing  factor  of  6.5  xl0~4  without  any  radial  dependence,  and 
extends  from  ^153  to  ~362  pc  from  the  center.  We  have  also 
assumed  the  “low  oxygen”  abundances  adopted  by  Alexander 
et  al.  (2000)  in  order  to  be  able  to  compare  our  results  with 


theirs.3  Because  the  grain  physics  has  been  updated  in  the  most 
recent  version  of  CLOUDY  (ver.  96),  we  have  included  the 
presence  of  grains  in  the  models,  using  “Orion-type”  grains.4  The 
inclusion  of  grains  also  allows  us  to  compute  the  thermal  dust 
continuum  emission  from  the  ionized  components  (see  §  4). 

The  inner  and  outer  radii  of  the  emission  regions  of  the  two 
components,  21,  109,  153,  and  362  pc,  correspond  to  angular 
distances  of  about  0726,  174,  1"9,  and  4"5,  respectively.  Table  2 
reports  the  predicted  line  fluxes  of  the  three  AGN  models,  A,  B, 
and  C,  together  with  the  observed  line  fluxes:  the  line  fluxes  are 

3  The  adopted  gas  phase  chemical  abundances  in  logarithmic  form  are  as 
follows:  H:  0.00,  He:  -1.00,  Li:  -8.69,  Be:  -10.58,  B:  -9.21,  C:  -3.43, 
N:  -3.96,  O:  -3.57,  F:  -7.52,  Ne:  -3.96,  Na:  -5.67,  Mg:  -4.43,  Al:  -5.53, 
Si:  -4.46,  P:  -6.49,  S:  -4.79,  Cl:  -6.72,  Ar:  -5.60,  K:  -6.88,  Ca:  -5.64, 
Sc:  -8.83,  Ti:  -6.98,  V:  -8.00,  Cr:  -6.33,  Mn:  -6.54,  Fe:  -4.40,  Co:  -7.08, 
Ni:  -5.75,  Cu:  -7.79,  and  Zn:  -7.40. 

4  The  abundances  of  the  grain  chemical  composition,  in  logarithmic  form,  are 
as  follows:  C:  -3.6259, 0:  -3.9526,  Mg:  -4.5547,  Si:  -4.5547,  andFe:  -4.5547. 


126 


SPINOGLIO  ET  AL. 


Fig.  1 . — AGN  ionizing  continua  used  as  input  for  the  photoionization  models 
of  NGC  1068.  The  three  continua  differ  in  the  frequency  region  in  the  range 
1  >  £Ryd  >  100,  while  outside  this  region  the  Pier  et  al.  (1994)  spectrum  was 
adopted.  The  solid  line  shows  the  continuum  derived  from  Alexander  et  al. 
(2000),  the  dashed  line  shows  a  simple  power-law  interpolation,  and  the  dotted 
line  shows  the  presence  of  the  predicted  big  blue  bump. 

given  for  each  of  the  two  components  1  and  2,  which  are  treated 
as  independent,  and  the  total  flux  for  each  model  is  simply  the 
sum  of  the  fluxes  of  the  two  components. 

We  can  see  from  Table  2  that  only  the  AGN  A  and  B  models, 
and  not  the  AGN  C  model,  reproduce  most  of  the  observed  high- 
ionization  line  fluxes.  The  low-  and  intermediate-ionization  lines 
are  expected  to  have  partial  or  frill  contributions  from  starburst 
and  PDR  components  (see  §  3.2).  This  first  result  rules  out  the 
presence  of  a  big  blue  bump  in  the  ionizing  continuum  of 
NGC  1068.  To  be  able  to  compare  the  modeled  ultraviolet  and 
optical  lines  with  the  observations,  we  also  listed  in  Tables  2  and 
3  their  dereddened  fluxes,  assuming  two  values  for  the  extinction: 
Eb^v  =  0.4  mag  (Malkan  &  Oke  1983)  and  EB_V  =  0.2  mag 
(Marconi  et  al.  1996).  We  find  that  the  AGN  B  model  over¬ 
predicts  several  of  the  intermediate-ionization  lines,  such  as  [S  iv] 
10.5  gm,  [Ne  in]  15.6  gm,  and  [S  in]  18.7  gm,  and  this  dis¬ 
crepancy  increases  when  adding  the  starburst  component  because 
these  lines  are  also  copiously  produced  by  that  component  (see 
the  next  section).  On  the  other  hand,  the  [Ne  n]  12.8  /im  emis¬ 
sion  is  underpredicted,  so  much  so  that  even  with  the  inclusion  of 
the  starburst  component  it  cannot  be  reproduced  with  this  model. 
As  we  discuss  further  in  §  3 .3,  a  composite  AGN  /starburst  model 
using  AGN  model  A  reproduces  the  [Ne  n]  12.8  gm  emission 
better  than  the  other  composite  models. 

3.2.  Modeling  the  Starburst  Ring 

NGC  1068  is  known  to  emit  strong  starburst  emission  from 
the  ringlike  structure  at  a  radial  distance  of  15"- 16"  from  the 
nucleus  (total  size  of  ~3  kpc),  traced  for  example  by  the  Br 7 
emission  (Davies  et  al.  1998).  Mid-IR  line  imaging  observations 
of  NGC  1068  have  been  published  by  Le  Floc’h  et  al.  (2001) 
based  on  ISOCAM  circular  variable  filter  (CVF)  observations. 
They  presented  an  image  of  the  7.7  gm  polycyclic  aromatic  hy¬ 
drocarbon  (PAH)  feature  that  shows  constant  surface  bright¬ 
ness  above  the  4th  contour  near  the  nucleus.  This  suggests  that 


star  formation  is  occuring  in  the  direction  of  the  nucleus  so  that 
nuclear  spectra  will  include  some  emission  from  star  forma¬ 
tion.  In  the  case  of  the  SWS  observations  that  we  are  modeling 
(reported  by  Lutz  et  al.  2000),  three  apertures  were  used  at 
different  wavelengths  with  the  two  largest  also  including  por¬ 
tions  of  the  brighter  starburst  ring  (see  Table  1).  To  estimate 
how  much  of  the  starburst  emission  is  contained  in  the  different 
apertures  used  in  the  observations,  we  have  used  a  continuum- 
subtracted  image  in  the  6.2  gm  feature  produced  by  C.  Dudley 
(2004,  private  communication)  using  the  same  ISOCAM  CVF 
data  set  examined  by  Le  Floc’h  et  al.  (2001).  The  6.2  gm  feature 
is  more  isolated  than  the  7.7  gm  feature,  which  is  blended  with 
the  8.7  gm  feature  and  the  silicate  absorption  feature,  but  the 
image  compares  well  with  the  published  7.7  gm  image  although 
we  have  zeroed  out  the  residuals  in  a  3x12  arcsec2  region 
centered  on  the  nucleus.  Based  on  this  image,  the  SWS  14x20, 
14x27,  and  20x33  arcsec2  slits  contain  13%,  23%,  and  46% 
of  the  6.2  gm  flux  contained  in  the  LWS  beam,  respectively, 
without  correction  for  the  neglected  region  of  poor  residuals 
(oriented  at  45°  to  our  synthetic  SWS  slits).  Since  PAH  features 
are  thought  to  be  a  good  tracer  of  PDRs  and  their  associated 
startbursts,  we  adopt  these  percentages  in  our  model  predictions 
of  SWS  line  strengths  in  the  starburst  models  presented  in  this 
section. 

In  fitting  our  starburst  models  to  the  observations,  we  have 
computed  the  line  fluxes  at  Earth  of  each  centrally  illuminated 
emitting  cloud  and  then  determined  the  number  of  clouds  needed 
to  best  fit  the  observed  line  fluxes.  We  have  chosen  the  starburst 
synthesis  modeling  program  Starburst99  (Leitherer  et  al.  1999) 
to  produce  input  ionizing  SEDs  for  the  CLOUDY  photo¬ 
ionization  code.  Models  were  followed  to  temperatures  down  to 
50  K  to  include  the  PDR  components.  We  compared  the  pre¬ 
dictions  of  an  instantaneous  star  formation  law  with  those  of  a 
continuous  star  formation  law.  For  both  types  of  models  we 
adopted  an  age  of  5  Myr,  a  Salpeter  initial  mass  function  (IMF; 
a  =  2.35),  a  lower  cutoff  mass  of  1  M@,  an  upper  cutoff  mass  of 
100  Mq,  solar  abundances  (Z  =  0.020),  and  nebular  emission 
included.  These  ionizing  SEDs  are  shown  in  Figure  2  with  a  total 
mass  of  M  =  106  M0  for  the  instantaneous  model  and  a  star 
formation  rate  of  1  Me  yr  for  the  continuous  model.  These 
particular  ionizing  continuum  shapes  were  selected  because  they 
are  consistent  with  the  Bry  equivalent  width  observed  by  Davies 
et  al.  (1 998)  in  the  starburst  ring.  We  have  estimated  that  the  Bry 
equivalent  width  in  each  of  the  individual  regions  of  the  map  of 
Davies  et  al.  (1998)  is  in  the  range  110-180  A.  According  to  the 
Leitherer  et  al.  (1999)  models  (see  their  Figs.  89  and  90),  for  a 
value  of  log  [JF(Bry,  A)]  >  2  only  instantaneous  models  with 
ages  less  than  ~6  x  1 06  yr  are  allowed. 

We  report  in  Table  3  the  line  fluxes  predicted  for  six  differ¬ 
ent  centrally  illuminated  starburst  models,  choosing  the  above 
instantaneous  star  formation  SED  as  the  input  ionizing  con¬ 
tinuum  and  using  CLOUDY  with  densities  of  nH  =  10,  100, 
and  1000  cnr3,  ionization  parameters  of  log  U  =  —2.5,  and 
—3.5,  and  an  inner  cloud  radius  of  50  pc.  As  a  function  of  the 
adopted  density,  we  then  determined  the  following  numbers 
of  emitting  clouds  needed  to  fit  the  observations:  33,000,  3300, 
and  330  clouds  for  the  three  values  of  the  density,  respectively. 
The  fluxes  reported  in  Table  3  are  therefore  the  total  starburst 
line  fluxes  at  Earth  and,  together  with  the  nuclear  line  fluxes  of 
Table  2,  can  be  compared  with  the  observations. 

We  have  also  mn  models  with  the  continuous  star  fomiation  law 
presented  above,  but  we  do  not  list  their  results  in  Table  3,  because 
the  differences  in  the  line  flux  predictions  compared  with  the  instan¬ 
taneous  models  are  insignificant,  compared  with  the  effects  of 


TABLE  2 

Comparison  of  Observed  Line  Fluxes  with  AGN  Model  Predictions 


Flux  (10  13  ergs  s  1  cm  2) 

AGN  A  Model0  AGN  B  Model1*  AGN  C  Model6 

Line  ID  and  Wavelength  -  -  - 


((im)  Observed /Da/Db  Comp.  1  Comp.  2  Comp.  1  Comp.  2  Comp.  1  Comp.  2 

O  vi  22(0.1032+  0.1037) .  37.4/4334/402  32.6  +  19.8  2.28  +  1.94  8.26  +  5.36  0.31  +  0.27  56.7  +  31.2  3.41  +2.10 

(Lya)„  20.1215 .  101.8/3562./602  179  619  149  500  239  638 

Niv]  20.1487 .  5.1/103./22.9  25.2  6.64  36.3  6.28  59.6  47.5 

(C  iv)„  20.1549 .  39.7/790/177  142  77.6  110.9  56.5  286  329 

He  ii  20.1640 .  21.4/426./95.5  112  85.9  108  81.8  346  345 

[Nev]  20.3426 .  15. 7/95738.7  97.7  23.3  68.2  7.36  271  115 

[Ne  in]  22(0.3869+  0.3968) .  19.2/97 743.2  37.1  +  11.2  41.5  +  12.5  102.8  +  31.0  102.  +  30.7  81.6  +  24.6  143.  +  43.0 

He  II  20.4686 .  6.1/27.6/13  15.0  12.3  14.7  11.8  46.5  48.4 

[O  hi]  22(0.4959  +  0.5007) .  25679647496  86.9  +  262  99.5  +  294  197.  +  593  185  +  557  148  +  446  241  +  726 

[Si  vi]  21.96 .  8. 0/9. 2/8. 6  11.2  1.53  10.8  1.86  23.3  3.74 

[Si  vii]  22.48 .  8.3  6.41  0.1  9.12  0.18  11.0  0.22 

[Si  ix]  22.584 .  3.0  0.49  ...  0.84  ...  0.32 

[Mg  viii]  23.028  .  11  2.97  ...  4.80  ...  4.29 

[Si  ix]  23.936 .  5.4  0.89  ...  1.56  ...  0.56 

[Mg  iv]  24.487 .  7.6  3.85  7.98  4.29  6.92  9.80  19.5 

[Arvi]  24.529 .  15  10.3  4.0  13.4  3.10  26.2  20.1 

[Mg  vn]  25.503  .  13  9.68  0.13  11.4  0.075  20.0  0.45 

[Mg  v]  25.610 .  18  10.7  8.63  7.85  3.84  34.6  29.5 

[Arn]  26.985 .  13  4.56  7.53  0.93  1.36  0.87  1.63 

[Na  hi]  27.318  .  5.8  0.56  0.59  1.27  1.51  0.67  1.21 

[Nevi]  27.652 .  110  153.7  11.9  105.4  3.2  393  62.4 

[Fe  vii]  27.815 .  3.0  1.83  1.63  1.65  0.71  6.34  6.30 

[Arv]  27.902 .  <12  2.14  1.83  2.62  1.92  3.72  6.52 

[Navi]  28.611 .  <16  1.20  0.16  0.95  ...  2.95  0.65 

[Arm]  +  [Mg  vii]  28.991 .  25  4.63  +  11.7  7.34  +  0.17  5.97  +  14.1  12.4  +  0.10  4.83  +  23.8  10.7  +  0.58 

[Fe  vii]  29.527 .  4.0  7.50  7.0  6.75  3.08  25.8  26.6 

[Siv]  210.510 .  58  38.4  60.5  85.3  126  75.8  187 

[Ne  n]  212.813 .  70  5.86  20.8  2.20  4.44  1.62  3.26 

[Arv]  213. 102 .  <16  2.23  2.73  2.75  2.88  3.93  9.57 

[Nev]  214.322 .  97  91.4  66.8  83.2  28.7  270  273 

[Ne  hi]  215.555 .  160  44.4  52.0  76  110  42.0  86.9 

[S  iii]  218.713 .  40  25.8  56.0  21.7  75.7  18.4  64.5 

[Nev]  224.317 .  70  40.8  52.6  35.8  22.1  116  215 

[O  iv]  225.890  .  190  24.1  71.9  28.3  69.7  56.4  240 

[S  iii]  233.481  .  55  8.00  29.8  5.57  41.1  4.84  38.4 

[Si  ii]  234.814 .  91  12.9  22.6  6.01  24.2  5.28  29.4 

[Ne  III]  236.013  .  18  3.32  4.4  5.74  9.40  3.20  7.46 

[O  iii]  251.814 .  110  9.04  33.9  16.3  63.2  9.06  48.5 

[N  III]  257.317 .  51  2.23  11.4  2.83  18.9  1.88  16.1 

[O  i  ]  263.184 .  156  3.84  1.8  1.57  2.26  1.24  3.10 

[O  iii]  288.356 .  110  1.22  8.2  2.25  15.5  1.23  12.6 

[N  ii]  2121.897 .  30  0.28  0.87  0.07  0.40  0.06  0.52 

[Oi]  2145.525  .  12  0.24  0.14  0.09  0.18  0.07  0.24 

[C  n]  2157.741 .  220  0.50  1.58  0.18  0.97  0.14  1.09 


a  Dereddened  line  flux,  assuming  Eg-y  =  0.4. 
b  Dereddened  line  flux,  assuming  Eg-y  =  0.2. 

°  AGN  A  parameters:  component  1:  log  U  =  —  1,  log  n  =  4.  internal  radius  ~21  pc,  external  radius  —  1 09  pc,  ionizing  spectrum  from  Alexander  et  al.  (2000): 
component  2:  log  U  =  —2,  log  n  —  3.3,  internal  radius  —  1 53  pc,  external  radius  —362  pc,  ionizing  spectrum  from  Alexander  et  al.  (2000). 
d  AGN  B  parameters:  same  as  AGN  A  models,  but  with  the  ionizing  spectrum  from  Pier  et  al.  (1994). 
e  AGN  C  parameters:  same  as  AGN  A  models,  but  with  the  ionizing  spectrum  that  includes  a  big  blue  bump  (see  text). 


TABLE  3 

Comparison  of  Observed  Line  Fluxes  with  the  Ring  Starburst  Model  Predictions 

Flux  (10  -13  ergs  s-1  cm-2) 

Line  ID  and  Wavelength  - 


(/zm)  Observed  SBR  Aa  SBR  Bb  SBR  C°  SBR  Dd  SBR  Ee  SBR  Ff 

O  vi  22(0. 1032  +  0.1037) .  37.4/4334/402 

(Ufa),,  20.1215 .  101.8/3562./602  1330  865  997  875  1010  917 

Niv]  20.1487 .  5.1/103./22.9  71.9  0.20  85.8  0.020  93.7  0.03 

(C  iv)„  20.1549 .  39.7/790/177  795  0.70  858  0.80  924  0.87 

He  ii  20.1640 .  21. 4/426795.5  693  50.2  696  50.5  696  50.8 

[Nev]  20.3426 .  15. 7/95738.7  58.1  ...  69.6  ...  73.9 

[Ne  in]  22(0.3869  +  0.3968) .  19.2/97743.2  314.  +  94.7  31.1  +  9.4  309.  +  93.1  31.5  +  9.5  316.  +  95.7  32.  +  9.7 

He  ii  20.4686 .  6.1/27.6/13  97.3  7.4  97.7  7.4  97.3  7.4 

[0  hi]  22(0.4959  +  0.5007) .  25679647496  1330  +  3990  25.1  +  75.6  1383  +  4158  26.0  +  78.2  1412  +  4257  26.8  +  80.5 

[Si  vi]  21.96 .  8.0/9.2/8.6 

[Si  vn]  22.48 .  8.3 

[Si  ix]  22.584 .  3.0 

[Mg  viii]  23.028 .  11 

[Si  ix]  23.936 .  5.4 

[Mg  iv]  24.487  .  7.6  3.30  ...  3.50  ...  3.53 

[Arvi]  24.529 .  15  6.20  ...  7.72  ...  8.02 

[Mg  vn]  25.503 .  13 

[Mg  v]  25.610  .  18  3.10  ...  3.63  ...  3.73 

[Arn]  26.985 .  13  5.90  3.70  4.29  3.70  4.12  3.73 

[Na  III]  27.318  .  5.8  0.60  0.10  0.54  0.11  0.54  0.11 

[Nevi]  27.652  .  110  16.9  ...  21.1  ...  21.8 

[Fe  vii]  27.815  .  3.0  0.40  ...  0.48  ...  0.50 

[ArV]  27.902 .  <12  7.50  ...  8.84  ...  9.11 

[Navi]  28.611 .  <16 

[Arm]  +  [Mg  vn]  28.991 .  25  46.5  7.10  43.2  7.13  42.2  7.13 

[Fe  vii]  29.527  .  4.0  1.70  ...  2.04  ...  2.11 

[Siv]  210.510 .  58  312  0.80  346  0.84  333  0.81 

[Ne  ii]  212.813  .  70  7.80  5.50  5.44  5.38  5.38  5.61 

[Ar  v]  213.102 .  <16  21.9  ...  25.5  ...  24.9 

[Nev]  214.322 .  97  271  ...  315  ...  314 

[Ne  iii]  215.555  .  160  295  50.5  276  50.2  270  49.8 

[S  III]  218.713 .  40  180  19.2  178  20.3  203  23.7 

[Nev]  224.317 .  70  303  ...  346  ...  294 

[0  iv]  225.890  .  190  1660  1.90  1762  2.00  1419  1.62 

[S  in]  233.481  .  55  785  87.5  673  81.2  350  46.2 

[Si  ii]  234.814 .  91  126  54.8  89.4  49.5  41.9  29.4 

[Ne  iii]  236.013  .  18  52.8  8.90  49.2  8.91  47.5  8.74 

[0  hi]  251.814 .  110  2980  81.5  3217  87.5  2524  71.9 

[N  III]  257.317 .  51  1110  40.0  921  35.0  372.9  14.7 

[Oi]  263. 184  .  156  165  130  113  127  116  124 

[0  iii]  288.356 .  110  4980  138  3729  109  977  29.0 

[N  Ii]  2121.897 .  30  80.0  37.3  47.8  29.7  12.3  8.78 

[Oi]  2145.525  .  12  16.2  12.9  11.0  12.5  10.1  11.3 

[C  ii]  2157.741 .  220  739  465  181  172  67.6  46.5 


a  SBR  A  parameters:  log  U  =  — 2.5 ,  log  n  =  1.0,  ionizing  spectrum  from  Starburst99  with  instantaneous  starformation  law,  M  =  106  M0,  IMF:  =  2.35  Mup  — 
100  Mq,  Mow  =  1  Mq,  nebular  emission  included,  Z  =  0.020,  age  of  5  Myr.  The  integration  was  stopped  at  a  temperature  of  50K,  the  adopted  abundances  are  those 
relative  to  H  n  regions  and  grain  emission  is  included.  The  adopted  number  of  clouds  is  33,000. 

b  SBR  B  parameters:  log  U  —  —3.5,  log  n  —  1.0,  all  other  parameters  as  for  SBR  A.  The  number  of  clouds  adopted  is  33,000. 

c  SBR  C  parameters:  log  U  —  —2.5,  log«  =  2.0,  all  other  parameters  as  for  SBR  A.  The  number  of  clouds  adopted  is  3300. 

d  SBR  D  parameters:  log  U  —  —3.5,  log  n  —  2.0,  all  other  parameters  as  for  SBR  A.  The  number  of  clouds  adopted  is  3300. 

e  SBR  E  parameters:  log  U  =  —2.5,  log  n  —  3.0,  all  other  parameters  as  for  SBR  A.  The  number  of  clouds  adopted  is  330. 
f  SBR  F  parameters:  log  U  —  —3.5,  log  n  —  3.0,  all  other  parameters  as  for  SBR  A.  The  number  of  clouds  adopted  is  330. 
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log  i/[Hz] 

Fig.  2. — Starburst  SEDs  used  as  input  for  the  photoionization  models  of 
NGC  1068.  The  two  continua  are  taken  from  Leitherer  et  al.  (1999)  and  rep¬ 
resent  a  continuous  starburst  model  ( solid  line)  and  an  instantaneous  model 
( dashed  line),  both  with  ages  of  5  Myr. 

density  and  ionization  parameter,  as  can  be  seen  from  Table  3. 
This  result  is  not  surprising  because  the  ionizing  continua  of 
the  two  starburst  models  are  quite  similar  in  shape  and  the  total 
number  of  clouds  is  a  free  parameter.  We  have  also  tried  contin¬ 
uous  starburst  models  with  much  longer  ages  (10,  20,  and  100  x 
105 6 *  yr),  but  because  the  shape  of  the  ionizing  continuum  again 
does  not  change  significantly,  the  resulting  emission-line  spec¬ 
trum  was  indistinguishable  from  that  derived  from  the  models  with 
an  age  of  5  x  106  yr. 

In  all  models  the  abundances  were  those  typical  of  H  n  re¬ 
gions,5  and  grains  of  Orion-type6  are  included.  The  integration 
was  allowed  to  run  until  the  temperature  of  the  gas  in  the  cloud 
cooled  to  T  =  50  K  in  order  to  include  the  PDRs  present  at  the 
interfaces  of  H  ii  regions  and  molecular  clouds. 

It  is  clear  from  Table  3  that  the  models  with  the  higher  ion¬ 
ization  parameter  (log  JJ  =  —2.5)  can  easily  be  ruled  out  be¬ 
cause  their  emission  in  many  intermediate-ionization  lines  is  far 
too  high  (see,  e.g.,  [O  iv]  26  pm,  [O  hi]  51,88  pm,  and  [N  ra] 
57  pm).  Among  the  models  with  the  lower  ionization  parameter 
(log  U  =  —3.5),  we  can  exclude  model  SBR  F,  with  density 
«h  =  1000  cm~3,  because  it  underestimates  many  far-lR  lines 
that  are  not  strongly  emitted  by  the  active  nucleus  (namely, 
[Si  n]  35  pm,  [N  in]  57  pm,  [O  hi]  88  /im,  [N  ii]  122  pn i,  and 
[C  ii]  158  /im),  while  the  low-density  model  (SBR  B,  with  «h  = 
10  cm-3)  overpredicts  the  [C  n]  158  /im  line  by  a  factor  of  2 
relative  to  the  other  far-IR  lines  and  does  not  reproduce  the  [O  in] 
doublet  ratio.  Finally,  the  intermediate-density  model  (SBR  D, 
with  «h  =100  cm~3)  gives  the  best  fit  to  the  observed  lines, 


5  The  adopted  gas  phase  chemical  abundances  in  logarithmic  form  are  as 
follows:  H:  0.00,  He:  -1.02,  Li:  -10.27,  Be:  -20.00,  B:  -10.05,  C:  -3.52, 
N:  -4.15,  O:  -3.40,  F:  -20.00,  Ne:  -4.22,  Na:  -6.52,  Mg:  -5.52,  Al:  -6.70, 
Si:  -5.40,  P:  -6.80,  S:  -5.00,  Cl:  -7.00,  Ar:  -5.52,  K:  -7.96,  Ca:  -7.70, 
Sc:  -20.00,  Ti:  -9.24,  V:  - 10.00,  Cr:  -8.00,  Mn:  -7.64,  Fe:  -5.52,  Co:  -20.00, 
Ni:  -7.00,  Cu:  -8.82,  and  Zn:  -7.70. 

6  The  abundances  of  the  grain  chemical  composition  in  logarithmic  form  are  as 

follows:  C:  -3.3249,  O:  -3.6516,  Mg:  -4.2537,  Si:  -4.2537,  andFe:  -4.2537. 


taking  into  account  that  the  AGN  component  must  be  added  to 
reproduce  the  total  flux  as  is  shown  in  §  3.3. 

We  estimate  the  average  PDR  parameters  using  the  models 
of  Kaufman  et  al.  (1999)  and  the  contour  plots  in  Luhman  et  al. 
(2003),  the  measured  [C  n]  158  and  [O  i]  145  fim  line  fluxes  (but 
not  the  [O  i]  63  /im  line  flux,  which  may  be  affected  by  absorp¬ 
tion  and/or  shocks),  and  the  far-IR  flux  integrated  over  the  LWS 
spectrum,  which  we  find  to  be  1.3xl0~8  ergs  cm~2  s-1. Here  we 
assume  that  the  [C  n]  line  emerges  predominantly  from  PDRs 
due  to  the  strong  starburst,  rather  than  the  diffuse  ionized  medium. 
With  this  assumption,  the  average  PDR  gas  density  and  UV  ra¬ 
diation  field  are  «h2  ~  1000  and  Go  ~  300,  respectively.  We 
note  that  if  we  instead  assume  that  the  [C  n]  line  flux  is  domi¬ 
nated  by  the  diffuse  ionized  medium,  using  the  correction  factor 
estimated  by  Malhotra  et  al.  (2001),  we  obtain  a  similar  gas 
density  «h2  ~  1500  but  a  significantly  higher  interstellar  radia¬ 
tion  field  Go  ~  1500.  For  both  cases,  the  parameters  derived 
are  in  the  range  of  those  of  the  normal  galaxies  in  the  Malhotra 
et  al.  (200 1 )  sample,  consistent  with  the  assumption  that  most  of 
the  far-IR  flux  originates  in  the  starburst  ring. 

3.3.  Adding  the  Two  Components 

Summing  the  line  intensities  of  each  one  of  the  two  com¬ 
ponents,  the  composite  spectrum  of  NGC  1068  can  be  derived 
and  compared  with  the  observed  one.  We  have  chosen  three 
combinations  to  compute  the  composite  models,  each  one  with 
a  different  AGN  model,  while  we  adopted  the  starbust  model 
with«H  =  100  cm~3  and  log  JJ  =  —3.5:  ( 1 )  composite  model  1 
(CM1),  with  the  AGN  ionizing  continuum  as  suggested  by 
Alexander  et  al.  (2000)  (model  AGN  A);  (2)  CM2,  with  the 
original  Pier  et  al.  (1994)  ionizing  spectrum  (model  AGN  B); 
(3)  CM3,  with  the  hypothetical  bump  (model  AGN  C).  The 
results  of  these  three  composite  models  are  given  in  Table  4, 
compared  to  the  observed  and  dereddened  values,  assuming  the 
two  choices  for  the  extinction  (see  §  3.1).  We  also  show  the 
results  of  the  three  composite  models  in  a  graphical  way  in 
Figure  3,  where  the  modeled  to  the  observed  flux  ratio  is  given 
for  each  line  for  the  case  of  an  extinction  of  EB_y  =  0.2  mag. 

A  simple  \2  test  of  the  three  models  resulted  in  a  reduced 
X2  of  1 1 .6,  17.1,  and  177  for  the  three  models  CM1,  CM2,  and 
CM3,  respectively,  for  an  extinction  of  EB_y  =  0.4  mag,  while 
these  values  become  23, 46,  and  325  for EB~v  =  0.2  mag.  Thus, 
of  the  models  explored,  CM1  with  EB  .y  =  0.4  mag  provides 
our  best  fit  to  the  observations.  We  note  that  model  CM1  re¬ 
produces  the  line  fluxes  to  within  a  factor  of  2.0  on  average, 
with  a  standard  deviation  of  1.3. 

4.  THERMAL  CONTINUUM  SPECTRUM 

To  model  the  total  mid-  and  far-infrared  thermal  dust  contin¬ 
uum  of  NGC  1068  we  have  used  different  computations  for  each 
component.  The  modeled  emission  is  shown  for  the  individual 
components  and  combined  in  Figure  4.  The  thermal  dust  emis¬ 
sion  from  the  AGN  NLRs  and  the  starburst  regions  in  the  ring 
have  been  computed  using  theCLOUDY  photoionization  mod¬ 
els  described  in  the  previous  sections.  Specifically,  the  UV  con¬ 
tinuum  reprocessed  by  the  dust  present  in  both  NLR  components 
1  and  2  of  model  AGN  A  has  been  diluted  by  the  same  covering 
factors  that  affect  the  line  emission  (c  =  0.45  and  0.29  for  com¬ 
ponents  1  and  2,  respectively).  Similarly,  we  computed  the  con¬ 
tinuum  emission  from  our  best  starburst  model  (SBR  D).  While 
we  find  that  the  continuum  produced  in  this  way  for  the  AGN 
is  consistent  with  the  observed  mid-infrared  energy  distribution 
(Fig.  4),  accounting  for  about  half  of  the  observed  emission  at 
mid-infrared  wavelengths,  the  emission  from  dust  associated  with 
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TABLE  4 

Comparison  of  Observed  Line  Fluxes  with  Composite  Model  Predictions 


Line  ID  and  Wavelength 

(m) 

Flux  (10  13  ergs  s 

1  cm”2) 

Observed/Da  /Db 

CM1C 

CM2d 

CM3e 

O  vi  22(0.1032  +  0.1037) . 

37.4/4334/402 

56.6 

14.2 

93.4 

(Lya)„  20.1215 . 

101.8/3562./602 

1673 

1524 

1752 

N  iv]  20.1487 . 

5.1/1 03./22.9 

31.8 

42.6 

107 

(C  iv)„  20.1549 . 

39.7/790/177 

220 

168 

616 

He  ii  20. 1640 . 

21. 4/426795.5 

248 

240 

741 

[Ne  v]  20.3426 . 

1 5.7/95738.7 

121 

75.6 

386 

[Ne  hi]  22(0.3869  +  0.3968) . 

19.2/97743.2 

143 

307 

333 

He  ii  20.4686 . 

6.1/27.6/13 

34.7 

33.9 

102.3 

[O  iii]  22(0.4959  +  0.5007) . 

25679647496 

847 

1636 

1665 

[Si  vi]  21.96 . 

8.0/9.2/8.6 

12.7 

12.7 

27.0 

[Si  vn]  22.48 . 

8.3 

6.5 

9.3 

11.2 

[Si  ix]  22.584 . 

3.0 

0.5 

0.8 

0.3 

[Mg  vra]  23.028 . 

11 

3.0 

4.8 

4.3 

[Si  ix]  23.936 . 

5.4 

0.9 

1.6 

0.6 

[Mg  iv]  24.487  . 

7.6 

11.8 

11.2 

29.3 

[Ar  vi]  24.529 . 

15 

14.3 

16.5 

46.3 

[Mg  vn]  25.503 . 

13 

9.8 

11.5 

20.5 

[Mg  v]  25.610 . 

18 

19.3 

11.7 

64.1 

[Ar  n]  26.985 . 

13 

15.8 

6.0 

6.2 

[Na  iii]  27.318 . 

5.8 

1.3 

2.8 

2.0 

[Ne  vi]  27.652 . 

110 

166 

109 

455 

[Fe  vn]  27.815 . 

3.0 

3.5 

2.4 

12.6 

[Ar  v]  27.902 . 

<12 

4.0 

4.5 

10.2 

[Na  vi]  28.611 . 

<16 

1.4 

0.9 

3.6 

[Arm]  +  [Mg  vn]  28.991 . 

25 

31 

39.7 

47.0 

[Fe  vn]  29.527  . 

4.0 

14.5 

9.8 

52.4 

[S  iv]  210.510 . 

58 

99.7 

212 

263 

[Ne  ii]  212.813 . 

70 

32.0 

12.0 

10.3 

[Ar  v]  213.102 . 

<16 

5.0 

5.6 

13.5 

[Ne  v]  214.322 . 

97 

158 

112 

543 

[Ne  iii]  215.555  . 

160 

147 

236 

179 

[S  iii]  218.713 . 

40 

102 

118 

103 

[Ne  v]  224.317 . 

70 

93.4 

57.9 

331 

[O  iv]  225.890 . 

190 

98 

100 

298 

[S  in]  233.481  . 

55 

119 

128 

124 

[Si  ii]  234.814 . 

91 

85 

79.7 

84.2 

[Ne  iii]  236.013  . 

18 

16.6 

24.0 

19.6 

[O  iii]  251.814 . 

110 

130 

167 

145 

[N  in]  257.317 . 

51 

48.6 

56.7 

53 

[O  i]  263.184 . 

156 

133 

131 

131 

[O  iii]  288.356 . 

110 

118 

127 

123 

[N  ii]  2121.897 . 

30 

30.8 

30 

30.3 

[O  i]  2145.525  . 

12 

12.9 

12.8 

12.8 

[C  ii]  2157.741 . 

220 

174 

173 

173 

a  Dereddened  line  flux,  assuming  Eg-v  =  0.4. 
b  Dereddened  line  flux,  assuming  Eg_v  =  0.2. 
°  CM1  =  AGN  A  +  SBR  D. 
d  CM2  =  AGN  B  +  SBR  D. 
e  CM3  =  AGN  C  +  SBR  D. 


the  starburst  ionized  and  photodissociated  regions,  although  sim¬ 
ilar  in  shape  to  the  observed  continuum,  produces  only  ~20%  of 
the  far-IR  continuum.  We  have  indeed  performed  a  search  in 
parameter  space  by  varying  the  age  of  the  starburst  [(4-6)  x 
106  yr],  the  gas  density  (from  10  to  1000  cm-3),  and  the  radius  of 
the  emitting  clouds  (from  25  to  100  pc),  but  no  starburst  model 
that  could  reproduce  simultaneously  the  observed  line  and  far- 
infrared  continuum  emission  was  found. 

If  the  CLOUDY  models  correctly  reflect  conditions  in  both 
the  ionized  and  photodissociated  gas,  these  results  may  imply 
that  the  bulk  of  the  starburst  far-infrared  continuum  arises  from 


dust  that  is  mixed  with  neutral  gas  not  directly  associated  with 
the  ionized  or  photodissociated  gas.  However,  because  photo¬ 
ionization  codes  such  as  CLOUDY  have  not  been  used  in  the 
past  to  model  the  dust  continuum,  and  because  we  may  not  have 
fully  searched  parameter  space,  we  are  hesitant  to  overinterpret 
these  results  until  more  detailed  comparison  with  galactic  ion¬ 
ized  and  photoionized  regions  are  carried  out.  We  have  there¬ 
fore  described  this  starburst  thermal  dust  component,  following 
Spinoglio  et  al.  (2002),  in  terms  of  a  graybody  function  with  a 
temperature  of  T  =  34  K,  and  a  colder  graybody  component 
at  T  =  20  K  (see  Fig.  4).  These  are  graybody  functions  with  a 
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Fig.  3. — Comparison  of  the  composite  models  with  the  observations  is 
shown  as  the  ratio  of  modeled  to  observed  flux  ratio  for  each  line,  with  the 
ionization  potential  in  the  x-axis.  The  assumed  reddening  is  E(B  —  V)  =  0.2. 
Panels  from  top  to  bottom  are  CM1,  CM2,  and  CM3,  respectively.  The  short- 
dashed  lines  represent  flux  ratios  within  a  factor  of  3  either  way. 


steep  (J3  =  2)  dust  emissivity  law.  Assuming  a  spherical  shell 
with  radius  of  1 .5  kpc  and  thickness  of  0.3  kpc,  the  inferred  aver¬ 
age  FE  molecular  density  associated  with  the  34  K  component  is 
3.7  cnT3.  The  total  mass  is  1.2x  109  and  2x  109  M0  for  the  34 
and  20  K  components,  respectively.  These  estimates  are  in  rea¬ 
sonable  agreement  with  the  ^4  x  109  M0  derived  for  the  molec¬ 
ular  ring  from  CO  emission  (Planesas  et  al.  1991). 

As  pointed  out  above,  the  dust  associated  with  the  ionized 
NLR  components  1  and  2  of  model  AGN  A  does  not  quite  ac¬ 
count  for  the  total  mid-infrared  continuum  (Fig.  4).  Therefore, 
we  have  assumed  that  the  missing  mid-infrared  arises  from 
the  “neutral-nuclear  component,”  which  has  been  observed  in 
a  variety  of  molecular  lines  (e.g.,  Tacconi  et  al.  1994;  Usero 
et  al.  2004).  For  this  neutral  component  we  have  assumed  a 
total  gas+dust  mass  of M  =  2x  107  M0  (Heifer  &  Blitz  1995) 
and  modeled  the  expected  continuum  using  a  nonlocal,  spher¬ 
ically  symmetric,  radiative  transfer  code  (Gonzalez-Alfonso  & 
Cemicharo  1997,  1999).  The  molecular  nuclear  emission  has 
been  resolved  into  a  circumnuclear  disk  or  ring  (Schinnerer 
et  al.  2000),  which  here  is  roughly  modeled  as  a  dusty  spherical 
envelope  with  inner  and  outer  radii  of  3  and  200  pc,  respec¬ 
tively.  We  assume  an  AGN  luminosity  of  3 .7  x  1010  Le  ,  which  is 
a  factor  of  ~3  lower  than  the  total  AGN  luminosity,  simulating 
that  most  of  the  AGN  radiation  escapes  through  the  poles  of 
the  molecular  disk  and/or  is  absorbed  in  the  NLR,  and  is  thus 
not  able  to  heat  the  molecular  gas  (e.g.,  Cameron  et  al.  1993). 
The  dust  envelope  is  divided  into  a  set  of  spherical  shells  in 
which  the  dust  temperature  is  computed  assuming  that  heating 
and  cooling  are  equal.  We  used  a  standard  silicate/amorphous 
carbon  mixture  with  optical  constants  given  by  Draine  (1985) 
and  Preibisch  et  al.  (1993).  The  density  profile  was  assumed  to 
be  oc  r  fl,  with  (3  regarded  as  a  free  parameter. 
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Fig.  4. — (a)  SED  of  NGC  1068  and  model  fit.  ISO- SWS  fluxes  are  taken 
from  Lutz  et  al.  (2000).  The  model  fit  is  the  composition  of  (1)  the  NLR  com¬ 
ponents  1  and  2  of  model  AGN  A,  (2)  the  emission  from  the  molecular  nu¬ 
clear  region,  (3)  the  34  K  starburst  emission,  and  (4)  the  cold  20  K  component. 
(b,  c )  Dust  temperature  vs.  the  radial  position  and  radial  continuum  opacity  vs. 
wavelength  for  the  nuclear  molecular  region. 


The  resulting  mid-infrared  emission,  obtained  with  (3  =  1 ,  is 
shown  in  Figure  4.  The  nuclear  molecular  component  has  an 
averaged  H2  density  at  the  inner  radius  of  («/( H2))  =  500  cm~3 
and  a  column  density  of/V(H2)  =  2x  1022  cm~2.  Once  the  emis¬ 
sion  from  this  neutral  component  is  summed  up  with  the  emis¬ 
sion  predicted  for  the  NLR  components  1  and  2,  a  good  fit  to 
the  observed  emission  for  X  >  9  /im  is  obtained.  For  a  <  9  /im, 
the  mid-infrared  emission  is  underestimated,  suggesting  the 
presence  of  a  hot  component,  probably  very  close  to  the  central 
AGN,  which  is  not  included  in  our  models. 

5.  THE  OH  LINES 
5.1.  General  Remarks 

In  NGC  1068,  we  detect  three  of  the  OH  rotational  lines,  all 
in  emission.  As  shown  in  the  energy  level  diagram  of  Figure  5, 
two  of  them  are  fundamental  lines,  connecting  the  ground-state 
2H3/23/2  level  with  the  2H3/25/2  (the  in-ladder  119  /im  line) 
and  the  2Hi/2l/2  level  (the  cross  ladder  79  /tm  line).  The  third 
line  is  the  lowest  transition  of  the  2ni/2  ladder:  the  163  /im  line 
between  the  J  =  3/2  and  J  =  1/2  levels.  The  detected  line 
fluxes  are  given  in  Table  1 .  The  fact  that  these  three  lines  are 


132 


SPINOGLIO  ET  AL. 


Vol.  623 


600 


iaO 

0) 


200 


0 


OH 

Hi/2  n3/2 
—  9/2 


Fig.  5. — Energy  level  diagram  of  OH.  Rotational  levels  with  energies  up 
to  600  K  are  shown;  the  three  lines  detected  in  NGC  1068  are  indicated  with 
arrows,  as  well  as  the  35  and  53  /im  lines  that  could  play  an  important  role  in 
the  radiative  excitation.  The  wavelengths  are  indicated  in  /im.  A-doubling  is 
ignored  because  the  A-doublets  are  not  resolved  with  the  ISO  grating  resolution. 


all  in  emission  is  in  striking  contrast  with  the  OH  lines  observed 
in  other  bright  infrared  galaxies,  such  as  Arp  220  (Fischer  et  al. 
1999;  Gonzalez-Alfonso  et  al.  2004),  Mrk  231  (Harvey  et  al. 
1999),  NGC  253  (Bradford  et  al.  1999),  and  M82  (Colbert 
et  al.  1999),  in  which  the  119  /jm  fundamental  is  in  absorption. 
The  79  //m  line  is  sometimes  seen  in  emission  and  sometimes  in 
absorption;  the  163  /im  line  is  always  seen  in  emission.  In  ad¬ 
dition  to  the  detections,  the  ISO- LWS  and  SWS  observations 
provide  upper  limits  on  fluxes  of  the  other  four  lines  that  arise 
between  the  lowest  six  rotational  levels.  The  LWS  spectra  in  the 
vicinity  of  the  detected  lines  (and  of  one  of  the  upper  limits  )  are 
shown  in  detail  in  Figure  6  ( histograms ).  In  this  section  we  dis¬ 
cuss  the  physical  conditions  necessary  to  excite  these  lines,  their 
probable  location  within  NGC  1068,  and  detailed  model  fits  of 
their  fluxes.  A  comparison  between  the  observed  and  modeled 
line  fluxes  is  given  in  Table  5  and  shown  in  Figure  6. 

5.2.  The  Excitation  Mechanism  of  the  OH  Lines 

The  unique  OH  emission-line  spectrum  of  NGC  1068  can 
provide  a  powerful  way  to  help  discriminate  between  the  prop¬ 
erties  of  the  molecular  clouds  in  NGC  1068  and  the  clouds  in 
other  galaxies  in  which  OH  has  been  observed.  Before  de¬ 
scribing  our  detailed  radiative  transfer  calculations,  it  is  instruc¬ 
tive  to  discuss  some  conclusions  that  are  model-independent. 
The  emission  in  the  OH  Tim  5/2-312  line  at  1 19  /mi  cannot  be 
explained  by  absorption  of  far-infrared  photons  followed  by 
cascade  down  to  the  upper  H1/2  5/2  level  of  the  transition. 
Rather,  we  argue  that  collisional  excitation  dominates.  Figure  5 
shows  the  energy  level  diagram  of  OH.  There  are  only  two  possi¬ 
ble  paths  to  excite  the  119  /tm  line  via  absorption  of  far- infrared 
photons:  via  the  35  /jiri  and/or  the  53  /im  ground-state  lines. 
Excitation  by  either  of  these  routes  has  other  observable  conse¬ 
quences.  In  the  case  of  simple  radiative  cascading,  the  Einstein 
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Fig.  6. — Comparison  of  the  observed  OH  lines  and  the  model  results.  As 
indicated  in  the  upper  panel,  the  upper  modeled  spectrum  ( solid  lines )  corre¬ 
sponds  to  the  model  for  the  nucleus,  the  middle  one  {dotted  lines )  corresponds 
to  the  starburst  modeled  as  a  whole,  and  the  lower  one  {dashed  lines)  corre¬ 
sponds  to  the  starburst  modeled  as  an  ensemble  of  individual  clouds  (see  text 
for  details). 


TABLE  5 

Comparison  of  Observed  Line  Fluxes  with  Model  Predictions 
for  the  Nuclear  Region 

Flux  (10~12  ergs  s-1  cm-2) 

Wavelength  - 


(/im)  Observed  Modeled  Notes 


34 .  <0.1  —0.43  Absorption 

53 .  <1.2  —0.66  Absorption 

79 .  1.1  1.13 

84 .  <1.2  0.06 

98 .  <1.2  0.17 

119 .  1.3  1.60 

163 .  0.38  0.35 
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A-coefficients  of  the  lines  involved  in  the  corresponding  cascades 
are  such  that  if  the  35  /tm  absorption  path  were  responsible  for 
the  observed  119  /tm  line  flux,  then  the  OH  ni/2  5/2 — 3/2  line 
at  98.7  /im  would  be  approximately  5  times  stronger  than  the 
119  /tm  line,  while  the  98.7  /tm  line  is  not  detected.  Hence,  this 
possibility  is  ruled  out.  Similarly,  if  absorption  in  the  53  /tm 
line  were  responsible  for  the  observed  119  /tm  line  flux,  then  the 
163  /tm  line  would  be  about  5  times  stronger  than  the  119  /tm 
line,  which  it  is  not.  We  can  therefore  conclude  from  the  con¬ 
straints  provided  by  the  other  far-infrared  OH  lines  that  the 
119  /im  emission  line  is  not  the  result  of  radiative  absorption 
and  cascading.  The  implication  is  that  OH  excitation  through 
collisions  is  more  important  in  NGC  1068  than  in  the  other  ob¬ 
served  galaxies  and  therefore  that  the  gas  responsible  for  the 
observed  emission  in  the  1 1 9  /tm  line  resides  predominantly 
in  relatively  dense  and  warm  environments  in  comparison  with 
these  other  sources. 

The  other  two  observed  emission  lines,  unlike  the  119  /tm 
line,  need  not  be  collisionally  dominated.  In  the  case  of  the 
n1/2  3/2-1/2  163  /im  line,  the  most  likely  excitation  mecha¬ 
nism  is  absorption  of  photons  emitted  by  dust  in  the  53  and 
35  /im  lines  followed  by  radiative  cascade.  The  upper  level  of 
this  transition  is  270  K  above  the  ground  state  (Fig.  5),  so  exci¬ 
tation  through  collisions  is  expected  to  be  ineffective  in  this  line. 
The  excitation  mechanism  of  the  ni/2-n3/2  1/2-3/2  79  /tm 
line  could  be  a  mixture  of  collisional  and  radiative  pumping. 
The  upper  level  of  this  transition  is  182  K  above  the  ground 
state,  so  that  a  warm  and  dense  region  could,  at  least  partially, 
excite  the  line  through  collisions.  Nevertheless,  the  line  could 
be  also  excited  through  the  same  infrared  pumping  mechanism 
that  results  in  the  observed  163  /tm  line  emission. 

In  conclusion,  the  119  /tm  line  is  collisionally  excited, 
whereas  absorption  of  photons  emitted  by  dust  in  the  53  and 
35  /im  lines  probably  dominates  the  excitation  of  the  163  /tm 
line.  The  79  /tm  OH  line  may  in  principle  be  excited  through 
both  mechanisms. 

5.3.  Constraints  on  the  Spatial  Origin  of  the  119  /tm  OH  Line 

In  NGC  1068  two  regions  with  very  different  physical  con¬ 
ditions  can  account  for  the  observed  OH  emission  as  discussed 
above:  the  compact  nuclear  region,  and  the  ring  and  bar  where 
the  starburst  is  taking  place.  A  warm  and  dense  region  is  re¬ 
quired  to  account  for  the  observed  1 1 9  /tm  line  emission,  given 
that  the  line  is  collisionally  excited,  so  the  warm  and  dense 
neutral  nuclear  region  around  the  AGN  should  be  considered  a 
good  candidate,  despite  its  small  size  (~5";  e.g.,  Planesas  et  al. 
1991;  Schinnerer  et  al.  2000),  for  the  following  reasons: 

1.  It  is  warm.  There  are  ~  1 0  ’  M.  of  hot  H2  (~2000  K) 
distributed  over  ~5"  (Blietz  et  al.  1994)  that  is  thought  to  be 
UV-  or  X-ray-heated  (Rotaciuc  et  al.  1991).  Both  PDRs  and 
X-ray-dominated  regions  (XDRs)  can  produce  a  range  of  tem¬ 
peratures  as  high  as  a  few  times  103  K  (Kaufman  et  al.  1999; 
Sternberg  &  Dalgamo  1995;  Maloney  et  al.  1996).  From  the 
CO  (4-3)  to  (1-0)  line  intensity  ratio,  Tacconi  et  al.  (1994)  de¬ 
rive  ^80  K  for  the  bulk  of  the  molecular  gas,  with  a  mass  of 
~3  x  107  M0  enclosed  in  within  the  central  4"  (Heifer  &  Blitz 
1995).  Lutz  et  al.  (2000)  have  reported  the  detection  of  pure 
H2  rotational  lines  within  the  ISO- SWS  aperture  and  estimated 
~2.5  x  107  M0  at  ^200  K,  but  these  lines  may  also  arise,  at  least 
partially,  from  the  inner  regions  of  the  3  kpc  starburst  ring. 

2.  The  molecular  clouds  within  the  nuclear  region  are 
dense,  although  there  is  some  dispersion  in  the  values  derived 
by  several  authors  based  on  HCN  emission:  Tacconi  et  al.  (1994) 


derived  an  H2  density  of  ~  1 05  cnr3,  whereas  subsequent  ob¬ 
servations  and  analysis  by  Heifer  &  Blitz  (1995)  yielded  a  den¬ 
sity  of  ~4  x  106  cm'3 .  An  intermediate  density  of  ~5  x  105  cm'3 
from  HCN  and  CS  (and  lower  for  other  tracers)  has  been  recently 
derived  by  Usero  et  al.  (2004). 

3.  The  OH  abundance  is  expected  to  reach  high  values  in 
regions  exposed  to  strong  incident  UV  fields  (PDRs;  Sternberg 
&  Dalgamo  1995)  and,  in  particular,  in  XDRs  (Lepp  &  Dalgamo 
1996).  The  remarkable  chemistry  found  by  Usero  et  al.  (2004) 
in  the  circumnuclear  disk  of  NGC  1068  is  indicative  of  an  over¬ 
all  XDR  and  suggests  a  high  OH  abundance  in  the  nuclear  region. 

Given  that  the  OH  119  /tm  line  is  collisionally  excited,  the 
possibility  that  the  line  might  arise  from  the  nuclear  region  can 
be  checked  by  computing  the  amount  of  warm  gas  required  to 
account  for  the  observed  emission: 


Mw{Mq)  =  1.6  xlO7 


10'5  ' 

5  x  105  cm  3 

Y(OH) 

[  «(H2)  J 

4.3  x  10  11  cm3  s  1 


(Clu) 


(1) 


where  Y(OH)  is  the  OH  abundance  relative  to  H2,  and  (c]u) 
is  the  collisional  rate  for  excitation  from  the  ground  Tim  3/2 
level  to  the  n2/2  5/2  one.  Equation  (1)  assumes  that,  although 
the  line  could  be  optically  thick,  it  is  effectively  optically  thin  and 
makes  use  of  the  observed  flux  of  1.2x1 O' 12  ergss'1  cm'2.  The 
reference  value  for  the  collisional  rate,  (ciu)  =  4.3  x  10'1 1  cm3 
s'1,  corresponds  to  gas  at  80  K  (Offer  et  al.  1994);  it  decreases 
by  a  factor  of  «2.7  for  gas  at  50  K  and  increases  by  a  factor  of 
3  for  gas  at  200  K. 

The  reference  OH  abundance  we  use  in  this  estimate,  1  O'5,  is 
the  result  of  two  separate  studies:  first,  calculations  of  molec¬ 
ular  abundances  by  Lepp  &  Dalgamo  (1996)  have  shown  that 
the  OH  abundance  in  XDRs  is  expected  to  be  about  2  orders  of 
magnitude  higher  than  the  abundance  of  HCN  and  HCO+.  The 
authors  in  fact  suggested  the  possibility  that  the  high  HCN/CO 
ratio  observed  in  the  nuclear  region  of  NGC  1068  could  be  a 
consequence  of  enhanced  X-ray  ionization.  Second,  the  pos¬ 
sibility  of  a  chemistry  dominated  by  X-rays  has  found  support 
from  observations  by  Usero  et  al.  (2004),  who  derive  abun¬ 
dance  ratios  of  HCN,  HCO+,  and  CN  in  general  agreement  with 
predictions  for  XDRs.  Since  the  HCN  abundance  derived  by 
Usero  et  al.  (2004)  is  ^10'7,  X(OH)  in  XDRs  could  attain 
values  as  high  as  10'5.  On  the  other  hand,  the  density  of  5  x 
105  cm'3  derived  by  Usero  etal.  (2004)  has  been  adopted  as  the 
reference  value  in  equation  (1).  Finally,  the  mass  of  molecular 
gas  derived  from  the  emission  of  several  molecular  tracers  is 
expected  to  be  at  least  ~2  x  1 07  M@  (Heifer  &  Blitz  1995),  which 
is  similar  to  the  value  required  in  equation  (1).  From  these  esti¬ 
mates  we  conclude  that,  if  the  OH  abundance  is  as  high  as  ~  1 0'5 
(i.e.,  if  the  predictions  for  XDRs  are  applicable  to  the  nuclear 
region  of  NGC  1068),  the  bulk  of  the  OH  119  /tm  line  could 
arise  there.  This  possibility  would  naturally  explain  why  NGC 
1068  is  unique  in  its  119  /tm  line  emission  among  galaxies  with 
full  LWS  spectra. 

Finally,  we  ask  whether  the  OH  119  /tm  line  could  arise  from 
an  even  more  compact  region,  i.e.,  from  a  toms  with  a  spatial 
scale  of  1  pc  surrounding  the  central  AGN.  According  to  typical 
parameters  given  by  Krolik  &  Lepp  (1989),  a  toms  is  expected 
to  be  hot  (~  1 03  K),  could  have  densities  of  1 07  cm'3 ,  and  there¬ 
fore  a  mass  of  ^105  M0.  Also,  the  OH  abundance  is  expected  to 
be  very  high,  5  x  10'5  to  5  x  10'4.  Equation  (1)  shows  that  the 


134 


SP1N0GLI0  ET  AL. 


Vol.  623 


relatively  low  mass  of  the  torus  (about  2  orders  of  magnitude 
lower  than  the  entire  nucleus)  could  be  compensated  by  the 
higher  density,  OH  abundance,  and  temperature  expected  there, 
so  this  possibility  cannot  be  neglected. 

The  reference  values  for  the  nuclear  abundance  and  density 
adopted  in  equation  (1)  are  rather  uncertain  (and  possibly  ex¬ 
treme).  The  continuum  models  of  §  4  indicate  that  the  mass  as¬ 
sociated  with  the  34  K  dust  component,  which  is  identified  with 
the  starburst  ring,  is  1.2x  109  M0.  If  ^5%  of  this  mass  corre¬ 
sponds  to  warm  molecular  gas  rich  in  OH,  the  amount  of  ex¬ 
tended  warm  gas  is  ~6  x  1 07  Me .  According  to  equation  ( 1 ),  the 
OH  emission  at  1 19  //m  can  then  also  be  explained  as  arising  in 
the  ring  if  the  associated  PDRs,  with  assumed  OH  abundance  of 
2  x  1 0  9  (Sternberg  &  Dalgamo  1995;  Goicoechea  &  Cemicharo 
2002;  Gonzalez-Alfonso  et  al.  2004),  have  densities  of  a  few 
times  105  cm~3.  Since  Papadopoulos  &  Seaquist  (1999)  found 
that  most  of  the  extended  molecular  gas  resides  in  dense,  com¬ 
pact  clouds,  this  scenario  also  seems  possible.  However,  the 
continuum  from  the  starburst  at  119  pm  is  strong,  so  one  ex¬ 
pects  that  equation  (1)  is  in  this  case  underestimating  Mw,  and 
the  quoted  physical  parameters,  X(OH )  and  n( H2),  are  lower 
limits.  The  effect  of  dust  emission  is  discussed  in  detail  below. 

In  conclusion,  a  definitive  answer  to  the  issue  of  the  spa¬ 
tial  origin  of  the  OH  119  pm  emission  cannot  be  inferred  from 
only  the  flux  observed  in  the  119  / tm  line.  Nevertheless,  useful 
constraints  on  this  subject  are  given:  the  line  could  either  be  ex¬ 
plained  as  arising  from  the  nucleus,  with  a  required  OH  abun¬ 
dance  ~  1 0  5 ,  or  from  the  extended  ring,  with  OH  abundance 
>2  x  1 0~6  and  density  greater  than  a  few  times  1 05  cnr3 .  Never¬ 
theless,  the  radiative  transfer  models  described  below,  which 
take  into  account  the  effect  of  the  continuum  emission  and  the 
excitation  of  the  79  and  163  /mi  lines,  point  toward  a  nuclear 
origin  of  the  OH  emission. 

5.4.  Outline  of  the  Models 

Analysis  of  the  OH  79  and  163  pm  lines  requires  the  use 
of  detailed  radiative  transfer  calculations  since,  as  pointed  out 
above,  the  emission  in  these  lines  is  expected  to  be  strongly  in¬ 
fluenced  by  absorption  of  far-infrared  continuum  photons.  We 
therefore  proceeded  to  model  the  OH  lines  with  two  different 
codes  and  confirmed  that  the  results  were  in  good  agreement 
with  each  other.  One  of  them,  described  in  Gonzalez-Alfonso 
&  Cemicharo  ( 1 997, 1 999),  has  recently  been  used  to  model  the 
far-infrared  spectrum  of  Arp  220  (Gonzalez-Alfonso  et  al. 
2004).  The  other  is  a  Monte  Carlo  radiative  transfer  code  used 
as  part  of  a  detailed  study  of  all  the  OH  lines  observed  by  ISO  in 
galaxies  (Smith  2004;  H.  A.  Smith  et  al.  2005,  in  preparation). 
The  code  was  developed  for  the  Submillimeter  Wave  Astronomy 
Satellite  (SWAS)  mission  and  is  a  modified  and  extended  ver¬ 
sion  of  the  Bemes  (1979)  code;  however,  it  includes  dust  as  well 
as  gas  in  the  radiative  transfer  and  also  corrects  some  optical 
depth  calculations  from  the  original  code  (Ashby  et  al.  2000). 
Both  methods  are  nonlocal,  non-LTE,  assume  spherical  sym¬ 
metry,  and  include  a  treatment  of  continuum  photons  from  dust 
mixed  in  with  the  gas.  Also,  both  codes  take  input  as  a  series  of 
concentric  shells,  each  of  which  is  assigned  a  size,  gas  and  dust 
temperature,  H2  density,  velocity  and  turbulent  velocity  width, 
and  molecular  abundance  relative  to  H2.  The  statistical  equi¬ 
librium  populations  of  OH  in  each  spherical  shell  are  computed 
by  including  the  excitation  by  dust  emission,  excitation  through 
collisions,  and  effects  of  line  trapping.  We  ran  two  models  to 
simulate  the  nucleus  of  the  galaxy  and  the  starburst  extended 
ring,  described  in  §§  5.4.1  and  5.4.2. 


5.4.1.  Models  for  the  Nuclear  Emission:  Constraints  on  the  Spatial 

Origin  of  the  79  and  163  /.tm  OH  Lines 

We  present  models  for  the  nuclear  OH  emission  that  im¬ 
plicitly  assume  that  the  119  pm  emission  line  arises  from  the 
nuclear  region :  X(OH)  =  10-5  is  adopted,  as  well  as  densities 
~5  x  105  cm~3  for  the  bulk  of  the  emitting  gas.  By  assuming  a 
pure  nuclear  origin  for  the  119  /.tm  line,  we  check  whether  the 
other  two  OH  lines  could,  in  such  a  case,  arise  from  the  same 
nuclear  region  or  require  a  more  extended  spatial  origin. 

The  models  use  the  dust  parameters  derived  for  the  nuclear 
molecular  region  described  in  §  4.  The  possible  contribution  to 
the  OH  excitation  of  far-infrared  photons  arising  from  the  com¬ 
ponents  1  and  2  of  model  AGN  A  is  ignored;  i.e.,  only  the  dust 
coexistent  with  the  molecular  gas  is  taken  into  account,  with  pre¬ 
dicted  flux  of  24-27  Jy  at  35-53  /.tm  (Fig.  4).  The  gas  temper¬ 
ature  is  assumed  to  be  uniform  and  equal  to  70  K  (e.g.,  Tacconi 
etal.  1994). 

The  H2  densities  derived  from  the  dust  model  [i.e.,  a  peak 
density  of  («/(H2))  =  500  cm~3  at  the  inner  radius]  are  not 
compatible  with  the  densities  inferred  from  different  molecular 
tracers.  This  indicates  that  the  medium  is  extremely  clumped,  as 
has  also  been  argued  elsewhere  (e.g.,  Cameron  et  al.  1993; 
Tacconi  et  al.  1994).  In  order  to  account  approximately  for  this 
dumpiness  in  our  models,  the  following  strategy  is  adopted:  we 
use  the  “real”  h(H2)  ~  5  x  105  cm~3  values  for  the  bulk  of  the 
gas  and  compute  the  volume  filling  factor  fv  =  («(H2))/m(H2), 
where  the  average  value  is  that  inferred  from  the  dust  model. 
The  expected  abundances  of  OH  and  the  dust  relative  to  H2, 
W(OH),  and  W(  dust)  are  then  multiplied  by  f,  so  that  the  right 
OH  and  dust  column  densities  are  used  in  the  calculations  to¬ 
gether  with  the  right  density  values.  The  same  density  profile  r~ 1 
that  was  used  in  the  dust  model  is  adopted,  so  that f,  is  uniform 
throughout  the  nuclear  region. 

The  modeled  fluxes  are  convolved  with  the  ISO- LWS  grat¬ 
ing  resolution  and  are  compared  with  the  data  in  Figure  6. 
Solid  lines  show  the  results  for  the  nuclear  model  that  assumes 
W(OH)  =  10~5  and  f,  =  2x  1(R5,  the  latter  value  implying  a 
density  in  the  outer  regions  (where  the  bulk  of  the  emission  is 
generated)  of  5  x10s  cm-3.  The  radial  OH  column  density  is 
A(OH)  =  2xl017  cm-2.  Besides  the  119  pm  line,  the  model 
nearly  reproduces  the  emission  in  the  79  and  163  pm  lines  and 
is  consistent  with  the  upper  limits  given  in  Table  5. 

We  also  checked  the  excitation  mechanism  of  the  other 
OH  lines  by  generating  an  additional  model  with  the  same  param¬ 
eters  as  above  except  for  the  continuum  emission,  which  is  now 
turned  off.  Therefore,  in  this  model,  the  lines  are  excited  ex¬ 
clusively  through  collisions  with  H2.  The  resulting  flux  of  the 
119  pm  line  remains  unchanged  when  the  dust  emission  is  ig¬ 
nored,  confirming  that  the  line  is  collisionally  excited.  On  the 
other  hand,  the  flux  densities  of  the  79  and  1 63  pm  lines  decrease 
in  the  “pure-collisional”  model  by  factors  of  2  and  6,  re¬ 
spectively,  showing  that  the  emission  in  these  lines  is  much  more 
affected  by  radiative  pumping.  We  conclude  that,  if  the  OH 
abundance  in  the  nucleus  were  high  enough  to  account  for  the 
collisionally  excited  119  pm  line,  the  observed  fluxes  in  the  79 
and  1 63  pm  lines  can  also  be  explained  as  arising  in  the  same 
nuclear  region. 

5.4.2.  Models  for  the  Starburst  Emission 

Two  simple  different  approaches  have  been  used  to  model 
the  OH  emission  from  the  starburst.  First,  we  have  roughly  mod¬ 
eled  the  whole  starburst  region  as  a  spherical  shell  with  external 
radius  of  1.5  kpc,  thickness  of  0.3  kpc,  and  average  H2  density 
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(«/( H2)}  =  3.7  cm-3,  so  the  corresponding  continuum  emission 
is  reproduced  with  7^  =  34  K  (§  4).  As  shown  in  §  5.3,  the  OH 
119  p m  emission  requires  densities  of  a  few  times  105  cm'3,  so 
we  have  assumed  a  volume  filling  factor  fv  =  7.5  xlO'6  and 
therefore  a  “real”  density  m(H2)  =  5  x  105  cm  3.  The  kinetic  tem¬ 
perature  is  assumed  to  be  7*.  =  100  K,  and  the  OH  abundance 
is  X(OH)  =  (2x10 ~6)fh,  where  f ,  =  0.05  is  the  assumed  frac¬ 
tion  of  warm  gas.  The  result  of  this  model  is  shown  in  Figure  6 
(upper  dotted  lines).  The  119  pm  line  is  reproduced,  but  the  flux 
densities  of  both  the  79  and  163  pm  lines  are  strongly  under¬ 
estimated.  The  reason  is  that  the  model  implicitly  assumes  that  the 
continuum  emission,  responsible  for  the  excitation  of  those  lines, 
arises  from  a  very  large  volume,  so  the  radiation  density  is  weak 
and  has  negligible  effect  on  the  line  excitation. 

Since  the  OH  emission  is  expected  to  arise  from  compact, 
discrete  PDRs  in  the  vicinity  of  O  or  early  B  stars,  where  the 
continuum  infrared  radiation  density  is  expected  to  be  stronger 
than  assumed  above,  we  have  also  tried  an  alternative  approach, 
which  consists  of  modeling  an  individual  “typical”  cloud  of 
the  starburst.  We  first  model  the  continuum  from  an  individual 
cloud  by  assuming  a  central  heating  source  and  computing  the 
dust  equilibrium  temperatures  at  each  radial  position  that  result 
from  the  balance  of  heating  and  cooling.  The  continuum  model 
is  adopted  if  (1)  leaving  aside  a  scaling  factor  (Nc,  the  number  of 
clouds  in  the  ensemble),  the  resulting  SED  is  similar  to  that  of 
the  far-infrared  emission  of  the  starburst  (i.e.,  the  34  K  com¬ 
ponent  found  in  §  4),  and  the  value  of  Nc  is  determined  by  re¬ 
quiring  that  the  absolute  continuum  flux  from  the  ensemble 
of  clouds  is  equal  to  that  observed  for  the  34  K  component; 
and  (2)  we  require  that  the  total  mass  of  the  ensemble  does  not 
exceed  the  mass  inferred  from  the  nonnuclear  region  (<4  x 
109  M0).  Once  the  continuum  is  fitted,  calculations  for  OH  are 
performed  by  assuming  Tk  =  Td,  and  Y(OH)  =  2x10  6. 

Several  models  with  various  density  profiles  were  found  to 
match  the  above  two  requirements.  The  common  characteristic 
of  all  of  them  is  the  relatively  high  column  density  of  the  in¬ 
dividual  clouds,  A(H2)  >  1023  cm-2,  which  is  a  consequence 
of  the  low  effective  dust  temperature  (34  K)  of  the  SED.  The 
results  of  the  simplest  model,  characterized  by  a  flat  density 
profile,  »( H2)  =  5xl05  cnr3,  are  given  here  for  reference. 
With  a  radius  of  5  x  1 0 1 7  cm,  a  stellar  luminosity  of  2  x  1 04  ZQ , 
andyVf  =  4  xlO6,  the  resulting  SED  is  similar  to  that  of  the  34  K 
component  (§  4).  For  these  clouds  we  obtain  a  total  mass  of 
1.8x  109  M0.  The  predicted  OH  emission/absorption  is  shown 
in  Figure  6  ( lower  dashed  lines).  In  spite  of  the  relatively  high 
density  and  temperatures  (28-250  K)  throughout  the  cloud,  the 
119  pm  line  is  predicted  to  be  too  weak,  and  the  79  ^m  line  is 
predicted  in  absorption.  We  have  found  this  result  quite  general: 
in  models  in  which  the  OH  abundance  is  high  enough  and  the 
radiation  density  becomes  strong  enough  to  pump  the  163  pm 
emission,  the  continuum  at  79  and  119  pm  is  absorbed  by  OH 
and  the  predicted  emission  in  the  corresponding  lines  is  reduced. 
Models  that  assume  a  density  profile  of  r  1  generally  predict  the 
119  pm  line  in  absorption.  In  some  models  in  which  the  OH 
abundance  was  allowed  to  vary  with  radial  position,  the  79  pm 
line  was  predicted  in  emission  but  by  far  too  weak  to  account 
for  the  observed  flux  density. 

In  conclusion,  no  starburst  model  is  found  to  satisfactorily 
reproduce  the  emission  observed  in  the  three  OH  lines.  If  the 
local  infrared  radiation  density  is  strong  enough  to  pump  the 
163  pm  line,  the  other  two  OH  lines  are  expected  to  be  weak  or 
in  absorption.  Furthermore,  the  high  density  assumed  for  the 
starburst  region  would  produce  a  relatively  high  HCN/CO  in¬ 


tensity  ratio,  which  is,  on  the  contrary,  ^0.01  in  the  spiral  arms 
(Heifer  &  Blitz  1995).  Finally,  the  PDR  models  described  in 
§  3.2  indicate  a  density  of  (1-1.5) xlO3  cm'3,  i.e.,  a  density 
much  lower  than  that  required  to  account  for  the  flux  density  of 
the  OH  119  pm  line.  Therefore(and  despite  the  simplicity  of 
our  models),  taken  together,  the  analysis  of  the  OH  lines  and 
the  derived  PDR  parameters  indicate  that  the  bulk  of  the  OH 
emission  arises  from  X-ray-dominated  nuclear  regions. 

6.  CONCLUSIONS 

The  main  results  of  this  article  can  be  summarized  as  follows: 

1.  The  complete  far-infrared  (50-200  pm)  spectrum  of 
NGC  1068  has  been  observed  for  the  first  time.  The  far-infrared 
ISO- LWS  spectrum  has  been  complemented  with  the  mid- 
infrared  data  of  7<S'0-SWS  and  with  shorter  wavelength  (UV, 
optical,  and  near-IR)  data  from  the  literature  to  assemble  a  com¬ 
posite  atomic  spectrum  as  complete  as  possible  with  the  aim  of 
modeling  the  different  line  emission  components  at  work.  This 
approach  has  been  especially  necessary  because  of  the  poor 
spatial  resolution  of  the  ISO  spectrometers,  which  were  not  able 
to  spatially  separate  the  emission  components.  The  lines  have 
been  interpreted  as  arising  from  two  physically  distinct  com¬ 
ponents:  the  AGN  component  and  a  starburst  component,  the 
first  one  nuclear  and  the  second  one  located  in  the  ring  at  a  radius 
of  15"- 16"  from  the  nucleus.  Both  components  are  character¬ 
ized  by  the  presence  of  dust  grains,  producing  strong  continuum 
emission  in  the  mid-  and  far-infrared.  The  density  and  ionization 
parameter  of  the  ~5  x  106  yr  old  starburst  are  found  to  be  «h  ~ 
100  cm'3  and  log  U  =  —3.5,  respectively.  Three  composite 
models  have  been  computed  with  different  AGN  components: 
the  first  one  has  the  ionizing  continuum  as  derived  from  Alexander 
et  al.  (2000),  showing  a  deep  trough  at  energies  of  a  few  rydbergs; 
the  second  has  the  monotonically  decreasing  ionizing  continuum 
given  by  Pier  et  al.  (1994),  and  the  third  has  a  big  blue  bump.  Two 
values  of  the  visual  extinction  ( Eb~v  =  0.2  and  0.4)  have  been 
adopted  to  correct  the  optical  and  ultraviolet  line  fluxes  for  the 
reddening.  The  agreement  between  the  composite  model  with  an 
AGN  ionizing  continuum  characterized  by  the  deep  trough  sug¬ 
gested  by  Alexander  et  al.  (2000)  is  very  satisfactory,  taking  into 
account  both  the  simplicity  of  the  photoionization  models  chosen 
to  avoid  dealing  with  too  many  free  parameters  and  the  large 
number  of  lines  that  originate  in  different  physical  regimes.  The 
agreement  between  the  observed  spectrum  and  what  is  predicted 
using  the  canonical  ionizing  continuum  is  slightly  poorer,  while 
the  presence  of  a  big  blue  bump  is  ruled  out. 

2.  The  50-200  pm  continuum  has  been  modeled  using  dif¬ 
ferent  components  arising  from  both  the  nucleus  and  the  star- 
burst  ring.  For  the  nucleus,  we  have  combined  the  dust  emission 
from  the  ionized  components  in  the  NLRs  modeled  by  CLOUDY 
with  the  neutral  component  reproduced  by  the  radiative  transfer 
code  used  for  the  OH  molecular  emission.  For  the  starburst  ring, 
our  CLOUDY  modeling  of  the  ionized+PDR  components  could 
not  reproduce  the  far-infrared  emission  while  simultaneously 
fitting  the  far-lR  lines.  Instead,  we  fitted  the  observed  continuum 
by  a  neutral  molecular  component,  reproduced  by  two  graybody 
components  at  temperatures  of  20  and  34  K,  assuming  a  steep 
(/ 3  =  2)  dust  emissivity  law. 

3.  The  unique  OH  emission  in  the  119  pm  line  cannot  be 
explained  in  terms  of  OH  excitation  through  absorption  of  35 
and  53  pm  photons  emitted  by  dust  but  rather  is  collisionally 
excited.  This  indicates  the  presence  of  a  warm  and  dense  region 
with  high  OH  abundance.  A  simple  excitation  analysis  yields  two 
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main  alternatives  for  the  “spatial  origin”  of  the  observed  119  /i  m 
line  emission:  ( 1 )  the  nuclear  region,  with  2  x  1 07  M@  of  warm  gas 
(80  K),  an  average  density  of  n{  H2)  =  5  x  105  cm  3,  and  an  OH 
abundance  of  ^10^5;  and  (2)  the  starburst  region,  if  ^5%  of  the 
associated  mass  (^6x  107  M@)  is  warm  (^100  K),  dense  (a  few 
times  105  cm~3),  and  rich  in  OH  \X(OH)  ~  2x  10  6], 

4.  Radiative  transfer  models  that  simulate  the  emission/ 
absorption  in  all  the  OH  lines  have  been  performed  for  both  the 
nuclear  and  the  starburst  region.  The  models  for  the  nucleus 
quantitatively  account  for  the  emission  in  the  three  OH  lines  if 
the  nuclear  physical  conditions  pointed  out  above  are  assumed. 
On  the  other  hand,  no  starburst  model  is  found  to  match  the  three 
OH  lines  simultaneously,  because  the  strong  far-infrared  con¬ 
tinuum  tends  to  produce  absorption,  or  to  weaken  the  emission, 
in  the  OH  119  and  79  //m  lines  (as  observed  in  other  galaxies). 
Therefore,  although  some  contribution  from  the  extended  star- 
burst  cannot  be  ruled  out,  our  models  indicate  that  the  bulk  of 
the  OH  emission  arises  in  the  nuclear  region.  The  high  nuclear 
OH  abundance  required  to  explain  the  emission  strongly  sug- 
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